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ABSTRACT 

We found significant correlations between the arrival directions of the highest energy 
photons (E 1 > 10 GeV) observed by EGRET and positions of the BL Lac type 
objects (BL Lacs). The observed correlations imply that not less than three per cent of 
extragalactic photons at these energies originate from BL Lacs. Some of the correlating 
BL Lacs have no counterparts in the EGRET source catalog, i.e. do not coincide 
with strong emitters of gamma-rays at lower energy. The study of correlating BL 
Lacs suggests that they may form a subset which is statistically different from the 
total BL Lac catalog; we argue that they are prominent candidates for TeV gamma- 
ray sources. Our results demonstrate that the analysis of positional correlations is a 
powerful approach indispensable in cases when low statistics limits or even prohibits 
the standard case-by-case identification. 

Key words: BL Lacertae objects: general - method: statistical - gamma-rays: theory. 



1 INTRODUCTION 



Identification of point sources observed at different wave- 
lengths is a standard problem in astronomy. This problem 
becomes very difficult at high energies, when the angular res- 
olution decreases and a mere positional coincidence becomes 
insufficient for the identif ication (an example i s the EGRET 
catalog of point sources dHartman et al.lll999l) more than a 
half of which have no optical/radio counterparts). The pho- 
ton flux also decreases with energy. Already at EGRET en- 
ergies the separation of point sources from the diffuse back- 
ground becomes a challenging problem by itself. At even 
higher energies, the flux becomes so low that it is more ap- 
propriate to think of the data as a collection of individual 
photons. In this situation, the standard method which con- 
sists in finding local excesses of flux and identifying these 
excesses with astrophysical objects becomes inadequate. 

In this paper we argue that in the case of low flux the 
problem of identification can be approached by an alterna- 
tive method based on the statistical analysis of correlations 
between individual photons and a given catalog of astro- 
physical objects. The advantage of this method is that it 
may give positive identification even in the case when the 
average number of photons from a source is much less than 
one. The price to pay is the statistical character of the ob- 
tained information: one may establish with certainty that 
the catalog contains actual sources, but may, in general, 
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not be able to identify them individually, nor estimate the 
luminosity of a given source. This method has been previ- 
ously applied to the c a se of ultra-h ig h energy cosmic r ays by 
i Tinvakov fc Tkachevl (1200 ll 12002^: lOorbunov et alJ l|20021 . 
12004) : lGorbunov fc Troitskvl J2004) . 

Here we apply this method of identification to the cat- 
alog of individual photons with energies E-i > 10 GeV pub- 
lished recently bv lThompson. Bertsch fc O'Neall J2004) . We 
demonstrate that the arrival directions of these photons cor- 
relate with positions of optically bright BL Lacertae type 
objects (BL Lacs), that is the BL Lacs are emitters of hard 
gamma-rays. We then analyse the emission and spectral 
properties of the correlating objects at different wavelengths 
and discuss possible implications of this analysis, in partic- 
ular for the prospects of TeV observations. 



2 CATALOGS OF CANDIDATE SOURCES 
AND GAMMA-RAYS 

Blasars, and the BL Lacs in particular, are the active galax- 
ies with jets pointing to the Earth. Many of them are known 
to be st rong gamma-ray em it ters ( se e, e.g.. [ Hartm an et alJ 
Jl999l) : Ivon Montignv et alJ (11993) : iMattox et alJ Jl997ft . 
They are believed to host powerful astrophysical acceler- 
ators of very energetic particles. 

Current understanding of physical propeties of BL Lacs 
is far from being complete. The standard approach suggests 
the two-bump broadband spectral energy distribution (SED; 
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see, e.g.. iPadovani fc Giommll lll995l) : lFossati et alJ il998^ 1. 

The lower-frequency peak, whose position varies from the 
infrared (low-energy-peaked BL Lacs, or LBL) to the X-ray 
(high-energy-peaked, or HBL) band, is often well-measured 
and is believed to be caused by the synchrotron radiation of 
relativistic electrons. The second bump is probably due to 
the synchrotron photons scattering off the same relativistic 
electrons (self-Compton, SC); its position should be corre- 
lated with that of the first bump and varies from MeV (LBL) 
to TeV (HBL) gamma-rays. Due to a scarcity of the gamma- 
ray data, the SC bumps are generally studied much worse 
than the synchrotron ones. 

While the two-bump SED is inherent to the electron- 
powered jets, strong gamma-ray emission from BL Lacs is 
expecte d also in other m odels, e.g. in the 'proton blazar' 
model jMannheiml H 9931 where it inevitably accompanies 
the acceleration of protons. Not surprisingly, many of blasars 
have been identified as EGRET sources <|Hartman et alJ 
ll999l:lvon Montignv et alJll995t iMattox et alJll997T) . There 
are no reasons for the gamma-ray spectra of blasars to have 
a cutoff at the EGRET energies; at least some of the objects 
are likely emitters of photons at higher frequencies. In a few 
cases, this has been confirmed by the TeV observations. 

The Table of BL Lac's (Table II) in the catalog of 
quasa rs and active galactic nuclei dVeron-Cettv fc Veronl 
2003) consists of objects with different spectral properties 
which are divided into three classes - confirmed ('BL'), high- 
polarization ('HP'), and probable/possible ('BL?') BL Lacs. 
This division is made a ccording to several criteria (so me of 
which are discussed bv IVeron-Cettv fc Veronl (J2000) ) and 
may reflect important differences in physical properties of 
the objects. In our correlation analysis we test these three 
subclasses separately. Each of them we divide in addition 
into optically bright (V < 18 mag) and dim (V ^ 18 mag) 
parts. In this way we obtain three subsamples of optically 
bright BL Lac's: 

(1) The set of all confirmed 'BL' type objects with the 
visual magnitude V < 18 mag. This set of BL Lac's contains 
178 objects. 

(2) The set of all confirmed 'HP' type objects with the 
visual magnitude V < 18 mag. This set of BL Lacs contains 
47 objects. 

(3) The rest of the objects liste d in the Table II of the 
catalog llVeron- Cettv fc Verorl2003h with the visual magni- 
tude V < 18 mag. This set consisting of bright unconfirmed 
BL Lac's contains 81 object. 

The catalog of EGRET gamma-rays of the high- 
est energies (E-, > 10 GeV) contains 1506 events 
jThompson. Bertsch fc O'Neal200i) . To suppress the back- 
ground of Galactic gamma-rays we make a cut on the Galac- 
tic latitude, b > 10°. This reduces the number of events 
down to 613. 



3 PROCEDURE 

Our analysis is based on the ca lculation of the angular corre - 
lation function as described bv lTinvakov fc Tkachevl feOOll) . 
The statistical significance of correlation is estimated by 
testing the hypothesis that the highest energy gamma-rays 
observed by EGRET and candidate sources are uncorrelated. 
This is done as follows. For a given set of sources and the 



angle 5, we count the number of pairs source - gamma ray 
separated by the angular distance less or equal to <5, thus 
obtaining the data count. We then replace the real data by 
a randomly generated Monte-Carlo set of gamma-rays and 
calculate the number of pairs in the same way, thus obtain- 
ing the Monte-Carlo count. We repeat the latter procedure 
many times calling successful those tries when the Monte- 
Carlo count equals or exceeds the data count. The number of 
successful tries divided by the total number of tries gives the 
probability P(S) that the excess in the data count occured 
by chance. The smaller is this probability, the stronger (more 
significant) is the correlation. Th e validity of this straightfor- 
ward approach does not depend llTinvakov fc Tkachevl2004T) 
on the completeness of the catalog of the candidate sources 
provided simulated sets of events correctly represent the de- 
tector exposure. 

The Monte-Carlo events are drawn from the pool of 
events generated according to the EGRET exposure map 
(available at ftp://cossc.gsfc.nasa.gov/compton/data/ 
egret/high_level/combined_data/). The latter depends on 
energy; we adopt the map relevant for the highest energy 
range 4 GeV < E 1 < 10 GeV. In our case the energy range 
is even higher; however, the corresponding exposure map 
is not available. We expect that this does not significantly 
influence our results. 

The significance of correlations is determined by 
the probability P(S) evaluated at the optimum value of 
S which can be obtained by Monte-Carlo si mulations 
jTinvakov fc TkachevteOOlllGorbunov et al.l2004h . This op- 
timum value is usually close to the detector angular resolu- 
tion. The EGRET dete ctor has been carefully calibrated by 
iThompson et all (Il993ft : its angular resolution depends on 
both the energy and the inclination angle. Averaged angular 
dispersion contains a narrow component and a wide-angle 
tail and can be fitted at a given energy by four Gaussians. 
The radius of a circle contai ning 67% of the gamm a-rays 
depends on energy as follows (IThompson et alJIlQQ^ . 1 

M£ 7 Ko.50°(H|^) a534 (1) 

which gives an estimate for the angle S. Because of the com- 
plexity of the EGRET angular resolution we do not fix S by 
the Monte-Carlo simulation. Instead, we follow an alterna- 
tive approach which consists in choosing the optimum bin 
size 8 from the data and correc ting the corresponding sig- 
nificance by the penalty factor llTinvakov fc Tkachevll200ll . 
120041: iFinlev fc Westerhofj l2004f) . We will see in the next 
section that this approach allows for a simple estimate of 
significance. 



4 RESULTS 

4.1 Positional correlations 

In FigQwe present the probabilities P(S) for all three sets 
of optically bright objects. Sets (1) and (2) exhibit strong 

1 The EGRET detector has been calibrated at energies E-y ^ 
10 GeV. In what follows we assume that Eq. is valid at least 
up to E-y ~ 30 GeV. 
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Figure 1. P(S) for the sets of bright (V < 18mag) 'BL', 'HP' 
and unconfirmed BL Lacs ('BL?') . 

correlations with the EGRET gamma-rays at separation an- 
gles compatible with Eq. Q , while correlations with the set 
of unconfirmed BL Lacs are absent. 

For the set (1), the minimum value of P(5) is P w 
1CF 11 , and there are 10 events which contribute to correla- 
tions in the minimum (0.37 events expected as background 
from random coincidences). For the set (2), the minimum 
value is P = 6.2 x 10 -9 with 7 events contributing to corre- 
lations (0.23 events expected from random coincidences). 2 

To obtain the significance of correlations one may take 
the lowest value of P{5) and m ultiply it by the penal t y fac- 
tor c al culated as described bv iTinvakov fc Tkachevl (2001, 
12004) : iFinlev fc Westerhofj J2004h . This factor is iust the 
number of statistically independent "attempts" to find the 
lowest probability. In the case at hand it can be replaced by 
a conservative upper bound: the penalty factor for variation 
of a quantity cannot exceed the number of steps, that is for 
S varying between 0° and 3° in steps of 0.05°, the penalty 
factor is ^ 60. Multiplying the minimum probabilities in 
Fig. by 60 clearly would not affect our conclusions. 

The minimum of P(S) for BLs is at 5^ n = 0.2°, while 
the minimum for HPs is at S^h-i = 0.35° . More compact clus- 
tering of gamma-rays around BLs as comapared to HPs can 
be explained by the difference in the typical energy of cor- 
relating gamma-rays (here and in what follows we count the 
event as correlating if it is at angular separation 8 ^ <5 m i n for 
the respective set). As can be seen from Table 1, energies of 
photons which correlate with BLs are systematically higher 
as compared to those associated with HPs. It follows from 
Eq. that the corresonding average angular resolutions 
are {5$ V> } = 0.34°, {5 { £ P) } = 0.46°, that could explain the 
observed hierarchy 5„?„ < 6„h. 

J mm mm 

Cumulative distributions of energies of gamma-rays 
correlating with BL and HP are shown in Fig. |5] The 
Kolmogorov-Smirnov (KS) test gives P=1.4% for the prob- 
ability that both sets of photons are drawn from one and the 
same distribution. The significance is not high, so the sys- 
tematic difference in energies could have occured by chance. 
Nevertheless it gives a hint that 'BL' and 'HP' objects may 

2 We report the probability calculated from the data count and 
average Monte-Carlo count assuming the Poisson distribution. 
This is a sufficiently good approximation at small angular scales. 
The direct calculation of probabilities below 10 is not feasible. 
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Figure 2. Cumulative distribution of energies of correlating 



gamma-rays for the sets of 'BL' and 'HP' type objects. 
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Table 1. BL Lac's 


from the 


samples (1), 


(2) ('BL' and 



'HP' in column 't') and correlating gamma-ray events. In 
the col umn 'Id', 'E' indica tes that the object is an EGRET 
source iHartman et alll999f) (positional identification of the cor- 
responding events_wlUi_the3EG_sources has been claimed by 
iThompson. Bertsch fc O'Neal <2004h 'l. 'G' - that it is a GeV 
source iLamb fc Macomb 1997) and 'T' — that it is a TeV source. 
The visual magnitudes are presented in column 'V, the red- 
shifts are presented in column 'z' (the question mark througout 
the table indicates that the value is unknown); the column 'F5' 
pres ents the radio-flux at 5 GHz in Jy (V, z and F5 are taken 
from lVeron-Cettv fc VerorJ bOOfl l. Note that 3EG J0433+2908 
and TXS 0506+056 correlate with two gamma-rays each, while 
PKS 2155-304 correlates with a triplet. 



have physically different spectra of gamma-rays. This issue 
cannot be elaborated further, in particular, because of the 
lack of knowledge of the distances to the objects in Tabled 
Note that PKS 0208-512 has a redshift z=l implying that 
the observed photon had two times higher energy at the 
source. 

Optically faint objects, V > 18 mag, of all three types 
(BL, HP and unconfirmed BL Lac's) do not correlate with 
gamma-rays. For these subsets P(S) > 10% in the range of 
S compatible with Eq. p). 
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4.2 Physical properties of correlating BL Lacs 

The important question is which observational and/or in- 
trinsic properties may distinguish efficient high-energy gam- 
ma-ray emitters. To systematically study this issue we carry 
out the KS test for compatibility of the distributions of 
correlating (without any cut on V) and of all objects of 
the same type (BL, HP) with respect to various parame- 
ters. Unlike correlations discussed in the previous section, 
this study may be affected b y incompleteness of the catalog 
(Veron-Cettv fc Veronl 12003ft . so the results of this section 
should be interpreted with care. 

Some of the correlating BL Lacs are emitters of mul- 
tiple photons, see Table 1. It is possible to treat such ob- 
jects in two ways. In the first approach multiple emitters 
are treated on equal footing with the rest of correlating BL 
Lacs. In the second approach multiple emitters are given 
weight which equals the number of photons observed from 
them. For the sake of completeness we present the results of 
both approaches. 

In Fig. E]we compare the cumulative distributions of vi- 
sual magnitudes, radio-fluxes at 5 GHz and X-ray fluxes 3 at 
1 keV for correlating and all objects of types 'BL' and 'HP'. 
Displayed curves correspond to the second approach when 
multiple emitters are weighted according to the number of 
observed photons. This allows to see positions of multiple 
emitters within relevant distributions. The results of the KS 
test for the first approach are given in parentheses. 

In Fig. [I] we compare the cumulitive distributions with 
respect to spectral indices. The radio-to-optic oro, radio- 
to- X-rays oirx and optic-to- X-rays aox indices are defined 
as cxab = \g(vAFA/vBF B )/\g{vA/vB), where Fa is the flu- 
ency at a frequency va- These parameters reflect intrinsic 
properties of the objects and are important for understand- 
ing the acceleration mechanism operating inside the sources. 
In particular, aox effectively measures the position of the 
synchrotron bump in the blazar's SED. 

Important physical information is contained also in the 
gamma-ray SEDs of the objects, which are available for the 
EGRET sources only (two objects of the BL type and three 
HPs) . The EGRET spectral indices of these five objects are 
unusually small, indicating higher fluency at higher energies. 
This is fully consistent with identification of them as the 
sources of E-, > 10 GeV photons. 

The sources of the 'HP' class have spectral indices 
aox > indicating that most of them are the HBLs with the 
synchrotron bump in X-rays. The IC bump is then expected 
at TeV energies, in consistency both with the EGRET spec- 
tral indices (for the 3EG sources) and with our evidence 
for energetic photons for the rest of the sample. Two of the 
sources have indeed been already detected in TeV. 

Much more interesting trends can be seen in the 'BL' 
sample. Most of these objects have significantly negative 
aox and are of the LBL type, which may be also seen 
directly from the strong dominance of the optical and ra- 



BL 



HP 




^ The catalog iVeron-Cettv fc Veronll2003ft contains visual mag- 
nitudes and radio-fluxes; X-ray fluxes are taken from the 
HEASARC database (http://heasarc.gsfc.nasa.govl or calcu- 
lated based on the count rates given there. Note that the emission 
and spectra of BL Lacs vary strongly with time, so this study 
should be considered as indicative only. 



Figure 3. Cumulative distributions of observational character- 
istics for correlating (thick lines) and all (thin lines) objects of 
the types 'BL' (left panels) and 'HP' (right panels). Upper panels 
show distributions of X-ray fluxes (10 _12 mW/m 2 ), middle panels 
correspond to radio flux (Jy) and bottom panels show apparent 
optical V-magnitudes. P corresponds to the Kolmogorov-Smirnov 
probability to obtain the distribution of correlating objects as a 
statistical fluctuation of the distribution of all objects. 



dio over X-ray emission (see Fig. suggesting the syn- 
chrotron bump in the infrared. This usually corresponds 
to the sub-GeV IC bump, but the possible detection of 
E 1 > 10 GeV photons from them (and, in both available 
cases, the EGRET spectral indices) indi cate rising spect ra 
at GeV energies (cf. the SED of ON 231 in lGhisellinil J2004l) i. 
Taken at face value, this fact has two immediate conse- 
quences: firstly, these objects are good candidates for the de- 
tection in TeV (up to now, only one of the eight has been de- 
tected); secondly, the two-bump electron-blazar model may 
not work, in its classical form, for their SEDs. The latter fact 
could b e explained, for i nstance, by a special geometry of the 
source l|Bednareklll99fti ) or by a more complicated gamma- 
ray-emission model. Clearly, the evidence is insufficient to 
claim that the eight sources considered here represent a new 
class (e.g., 'proton blazars'). However, one may note that 
it is the 'BL' objects brighter than 18 mag which correlate 
with ultra- high-energy cosmic rays 4 , and the acceleration of 
protons should take place in some of them, inevitably ac- 
companied by the emission of energetic gamma-rays. This 
tempting conjecture awaits further studies with the existing 
(EGRET) and upcoming (GLAST) gamma-ray data. 

Another interesting feature is related to the doublets 



4 This statement relates to the whole sample of 178 objects. Of 
eight correlating BLs in Tab. only one falls in the errorbox of 
the arrival directions of cosmic rays - 3EG J0433+2908 coincides 
with an AGASA-detected doublet. 
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Figure 4. Cumulative distributions of intrinsic spectral prop- 
erties for correlating (thick lines) and all (thin lines) objects of 
the types 'BL' (left panels) and 'HP' (right panels). Upper, mid- 
dle and bottom panels show distributions of ctRx > a RO an d Q OX 
correspondingly. P corresponds to the Kolmogorov-Smirnov prob- 
ability to obtain the distribution of correlating objects as a sta- 
tistical fluctuation of the distribution of all objects. 



and the triplet of gamma-rays correlating with BL Lac's. 
All of them are monochromatic, see Table 1. This may hint 
at the details of the accelarating mechanism in the sources 
or ambient matter density and magnetic field strength close 
to the source. 

As we have already pointed out, the small number of 
corre lating photons and incomple teness of the original cat- 
alog dVeron-Cettv fc Veronll2003l) prevent one from making 
definite conclusions on the basis of the analysis performed 
in this subsection. However, our results suggest that the BL 
Lacs which are bright in radio and optical band may of- 
ten emit energetic gamma-rays. Together with the fact that 
correlating BL's on average emit more energetic photons (in 
E~t > 10 GeV band, see Fig. this may indicate a devia- 
tion from the conventional two-bump SED model, possibly 
related to the acceleration of cosmic rays. 



5 CONCLUSIONS 

The statistical analysis we performed in this paper reveals 
that the arrival directions of the high-energy EGRET pho- 
tons coincide with positions of BL Lacs ('BL' and 'HP' ob- 
jects) far too often to be explained by chance: out of 10 
coincidences observed for BLs only 0.37 (in average) would 
be expected if the photons and BLs were uncorrelated. Thus, 
our analysis has established with certainty that the BL cat- 
alog contains sources of gamma-rays. Moreover, nearly all 
correlating photons are due to sources, and therefore all cor- 
relating BLs are likely to be actual emitters (the same con- 



cerns HP objects). On the other hand, one may check that 
there is no significant clustering in the set of highest-energy 
EGRET photons, so the standard methods of identification 
would be inconclusive. 

The fact of correlations is not the only conclusion which 
statistical methods are able to establish. As has been shown 
in sect. 4.2, they allow to address a question of which phys- 
ical property characterizes the actual emitters. This, how- 
ever, requires the completeness of the catalog of candidate 
sources, as well as better statistics. 
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